function [omega, peaks] = find_modes(fourier_series, sam_len)

    [pk0,loc0] = findpeaks(fourier_series);
    [pk1,loc1] = findpeaks(pk0);
    freq = loc0(loc1);
    omega = freq * (2 * pi)/sam_len;
    peaks = pk1;

    return
    
end